VTrack has some cool tools including COA and estimating activity spaces from passive telemetry data. This package has a suite of functions that help users follow a workflow to efficiently setup, analyse and visualise telemetry data. For more details on the workflow check out this vignette.

Setup data

The primary object VTrack uses is called an ATT object. This object houses the three core datasets required to analyse telemetry data:

  1. Detections
  2. Receiver deployment metadata
  3. Tag metadata

Apart from these datasets, the ATT object also stores information on the Geographic Coordinate Reference System for that dataset. GLATOS is capable of exporting OTN data to a format readable by VTrack. Conversion of OTN data requires data from the OTN ERDDAP server Data can be found at the following links: Animals, Receivers, and Tag Releases.

Lets load some data from the glatos package, remove false detections and set it up to use in VTrack

library(tidyverse)
library(glatos)
library(plotly)
wal_det_file <- system.file("extdata", "walleye_detections.csv", package = "glatos")
walleye_detections <- read_glatos_detections(wal_det_file)

rec_file <- system.file("extdata", "sample_receivers.csv",  package = "glatos")
rcv <- read_glatos_receivers(rec_file)

# find and remove false detections
filtered_detections <-
  walleye_detections %>% 
  false_detections(tf = 3600) %>% 
  filter(passed_filter != FALSE)
## The filter identified 93 (1.3%) of 7180 detections as potentially false.
# Now we can convert the GLATOS filtered data into a ATT object to be used in VTrack
ATTdata <- glatos::convert_glatos_to_att(walleye_detections, rcv)
## [1] "http://www.marinespecies.org/rest/AphiaRecordsByVernacular/walleye"
attr(ATTdata, "CRS") <- sp::CRS("+init=epsg:4326")

head(ATTdata)
## $Tag.Detections
## # A tibble: 7,083 x 8
##    Date.Time           Transmitter Station.Name Receiver Latitude Longitude
##    <dttm>              <fct>       <fct>        <fct>       <dbl>     <dbl>
##  1 2012-04-29 01:48:37 A69-9001-3… TTB-002      VR2W-11…     43.4     -84.0
##  2 2012-04-29 01:52:55 A69-9001-3… TTB-002      VR2W-11…     43.4     -84.0
##  3 2012-04-29 01:55:12 A69-9001-3… TTB-002      VR2W-11…     43.4     -84.0
##  4 2012-04-29 01:56:42 A69-9001-3… TTB-002      VR2W-11…     43.4     -84.0
##  5 2012-04-29 01:58:37 A69-9001-3… TTB-002      VR2W-11…     43.4     -84.0
##  6 2012-04-29 02:01:22 A69-9001-3… TTB-002      VR2W-11…     43.4     -84.0
##  7 2012-04-29 02:03:47 A69-9001-3… TTB-002      VR2W-11…     43.4     -84.0
##  8 2012-04-29 02:05:33 A69-9001-3… TTB-002      VR2W-11…     43.4     -84.0
##  9 2012-04-29 02:08:00 A69-9001-3… TTB-001      VR2W-11…     43.4     -84.0
## 10 2012-04-29 02:08:00 A69-9001-3… TTB-002      VR2W-11…     43.4     -84.0
## # … with 7,073 more rows, and 2 more variables: Sensor.Value <int>,
## #   Sensor.Unit <fct>
## 
## $Tag.Metadata
## # A tibble: 3 x 12
##   Tag.ID Transmitter Common.Name Sci.Name Tag.Project Release.Latitude
##    <int> <fct>       <fct>       <fct>    <fct>                  <dbl>
## 1    153 A69-9001-3… walleye     Sander … HECWL                     NA
## 2     22 A69-9002-1… walleye     Sander … HECWL                     NA
## 3     23 A69-9002-1… walleye     Sander … HECWL                     NA
## # … with 6 more variables: Release.Longitude <dbl>, Release.Date <date>,
## #   Sex <fct>, Tag.Life <int>, Tag.Status <fct>, Bio <fct>
## 
## $Station.Information
## # A tibble: 898 x 9
##    Station.Name Receiver Installation Receiver.Project Deployment.Date    
##    <fct>        <fct>    <fct>        <fct>            <dttm>             
##  1 WHT-009      VR2W-10… <NA>         HECWL            2010-09-22 18:05:00
##  2 FDT-001      VR3-442  <NA>         DRMLT            2010-11-12 15:07:00
##  3 FDT-004      VR3-441  <NA>         DRMLT            2010-11-12 15:36:00
##  4 FDT-003      VR3-444  <NA>         DRMLT            2010-11-12 15:56:00
##  5 FDT-002      VR3-447  <NA>         DRMLT            2010-11-12 16:26:00
##  6 DTR-001      VR3-439  <NA>         DRMLT            2010-11-12 19:43:00
##  7 DTR-002      VR3-445  <NA>         DRMLT            2010-11-12 20:01:00
##  8 DTR-003      VR3-440  <NA>         DRMLT            2010-11-12 20:21:00
##  9 DTR-004      VR3-446  <NA>         DRMLT            2010-11-12 20:48:00
## 10 LVD-001      VR2W-10… <NA>         HECWL            2011-04-21 14:51:00
## # … with 888 more rows, and 4 more variables: Recovery.Date <dttm>,
## #   Station.Latitude <dbl>, Station.Longitude <dbl>, Receiver.Status <fct>
library(VTrack)
## Loading required package: foreach
## 
## Attaching package: 'foreach'
## The following objects are masked from 'package:purrr':
## 
##     accumulate, when
## Registered S3 methods overwritten by 'adehabitatMA':
##   method                       from
##   print.SpatialPixelsDataFrame sp  
##   print.SpatialPixels          sp
## Registered S3 method overwritten by 'xts':
##   method     from
##   as.zoo.xts zoo
aba <- abacusPlot(ATTdata, det.col = "red")

# Abacus plot by station
aba_facet <- abacusPlot(ATTdata, det.col = "red", facet = T)
print(aba)

print(aba_facet)

library(lubridate, verbose = FALSE, warn.conflicts = FALSE)

abacus_data <- 
  ATTdata$Tag.Detections %>% 
  left_join(ATTdata$Tag.Metadata,  by = "Transmitter") %>% 
  mutate(date = date(Date.Time)) %>% 
  group_by(Tag.ID, date, Sci.Name) %>% 
  summarise(Num.Det = n())

abacus_data %>% 
  ggplot(aes(x = date, y = factor(Tag.ID), size = Num.Det, color = Sci.Name)) +
  geom_point() +
  theme_linedraw() +
  labs(title = "Daily detection plot", x = "Date", y = "Tag ID", color = "Species", size = "Number of\ndetections")

detSum <- detectionSummary(ATTdata, sub = "%Y-%W")

detSum
## $Overall
## # A tibble: 3 x 10
##   Tag.ID Transmitter Sci.Name Sex   Bio   Number.of.Detec… Number.of.Stati…
##    <int> <fct>       <fct>    <fct> <fct>            <int>            <int>
## 1     22 A69-9002-1… Sander … FEMA… <NA>              2807               11
## 2     23 A69-9002-1… Sander … FEMA… <NA>              1327               23
## 3    153 A69-9001-3… Sander … FEMA… <NA>              2949               54
## # … with 3 more variables: Days.Detected <int>, Days.at.Liberty <dbl>,
## #   Detection.Index <dbl>
## 
## $Subsetted
## # A tibble: 35 x 12
##    Tag.ID subset Transmitter Sci.Name Sex   Bio   Number.of.Detec…
##     <int> <fct>  <fct>       <fct>    <fct> <fct>            <int>
##  1     22 2012-… A69-9002-1… Sander … FEMA… <NA>               205
##  2     22 2012-… A69-9002-1… Sander … FEMA… <NA>               453
##  3     22 2012-… A69-9002-1… Sander … FEMA… <NA>               163
##  4     22 2013-… A69-9002-1… Sander … FEMA… <NA>                30
##  5     22 2013-… A69-9002-1… Sander … FEMA… <NA>              1095
##  6     22 2013-… A69-9002-1… Sander … FEMA… <NA>               762
##  7     22 2013-… A69-9002-1… Sander … FEMA… <NA>                79
##  8     22 2013-… A69-9002-1… Sander … FEMA… <NA>                20
##  9     23 2012-… A69-9002-1… Sander … FEMA… <NA>               272
## 10     23 2012-… A69-9002-1… Sander … FEMA… <NA>               501
## # … with 25 more rows, and 5 more variables: Number.of.Stations <int>,
## #   Days.Detected <int>, New.Stations <int>, Days.at.Liberty <dbl>,
## #   Detection.Index <dbl>
# Overall metrics
head(detSum$Overall)
## # A tibble: 3 x 10
##   Tag.ID Transmitter Sci.Name Sex   Bio   Number.of.Detec… Number.of.Stati…
##    <int> <fct>       <fct>    <fct> <fct>            <int>            <int>
## 1     22 A69-9002-1… Sander … FEMA… <NA>              2807               11
## 2     23 A69-9002-1… Sander … FEMA… <NA>              1327               23
## 3    153 A69-9001-3… Sander … FEMA… <NA>              2949               54
## # … with 3 more variables: Days.Detected <int>, Days.at.Liberty <dbl>,
## #   Detection.Index <dbl>
# Weeky metrics
head(detSum$Subsetted)
## # A tibble: 6 x 12
##   Tag.ID subset Transmitter Sci.Name Sex   Bio   Number.of.Detec…
##    <int> <fct>  <fct>       <fct>    <fct> <fct>            <int>
## 1     22 2012-… A69-9002-1… Sander … FEMA… <NA>               205
## 2     22 2012-… A69-9002-1… Sander … FEMA… <NA>               453
## 3     22 2012-… A69-9002-1… Sander … FEMA… <NA>               163
## 4     22 2013-… A69-9002-1… Sander … FEMA… <NA>                30
## 5     22 2013-… A69-9002-1… Sander … FEMA… <NA>              1095
## 6     22 2013-… A69-9002-1… Sander … FEMA… <NA>               762
## # … with 5 more variables: Number.of.Stations <int>, Days.Detected <int>,
## #   New.Stations <int>, Days.at.Liberty <dbl>, Detection.Index <dbl>
detSum$Subset %>% 
  mutate(Date = as.Date(paste(subset, "01", sep = "-"), "%Y-%W-%w")) %>% 
  ggplot(aes(x = Date, y = Number.of.Detections, colour = factor(Tag.ID))) +
  geom_line() + geom_point() +
  labs(title = "Number of Detections over time", x = "Date", y = "Number of Detections", colour = "Tag ID") +
  theme_linedraw()

dispSum <- dispersalSummary(ATTdata)

head(dispSum)
## # A tibble: 6 x 24
##   Date.Time           Transmitter Station.Name Receiver Latitude Longitude
##   <dttm>              <fct>       <fct>        <fct>       <dbl>     <dbl>
## 1 2012-03-27 13:05:27 A69-9002-1… MAU-001      VR2W-11…     41.6     -83.6
## 2 2012-03-27 13:16:29 A69-9002-1… MAU-001      VR2W-11…     41.6     -83.6
## 3 2012-03-27 13:20:49 A69-9002-1… MAU-001      VR2W-11…     41.6     -83.6
## 4 2012-03-27 13:25:29 A69-9002-1… MAU-001      VR2W-11…     41.6     -83.6
## 5 2012-03-27 13:28:07 A69-9002-1… MAU-001      VR2W-11…     41.6     -83.6
## 6 2012-03-27 13:41:56 A69-9002-1… MAU-001      VR2W-11…     41.6     -83.6
## # … with 18 more variables: Sensor.Value <int>, Sensor.Unit <fct>,
## #   Tag.ID <int>, Common.Name <fct>, Sci.Name <fct>, Tag.Project <fct>,
## #   Release.Latitude <dbl>, Release.Longitude <dbl>, Release.Date <date>,
## #   Sex <fct>, Tag.Life <int>, Tag.Status <fct>, Bio <fct>,
## #   Release.Dispersal <lgl>, Release.Bearing <lgl>,
## #   Consecutive.Dispersal <dbl>, Consecutive.Bearing <dbl>,
## #   Time.Since.Last.Detection <dbl>
dispSum %>%
  mutate(Date = date(Date.Time),
         Tag.ID = factor(Tag.ID)) %>% 
  group_by(Date, Tag.ID) %>% 
  summarise(Max.Disp.km = max(Consecutive.Dispersal)) %>%
  ggplot(aes(x = Max.Disp.km, y=..scaled..,  fill = Tag.ID, color = Tag.ID)) +
  geom_density(alpha=0.1) +
  geom_rug(sides = "t", aes(y = Max.Disp.km, color = Tag.ID)) +
  ylim(0,1) +
  labs(x = "Maximum distance traveled per day (km)",
       y = "Relative frequency of daily movements",
       fill = "Tag ID", color = "Tag ID") +
  theme_linedraw()
## Warning: Removed 3 rows containing non-finite values (stat_density).

library(sf)
## Linking to GEOS 3.7.2, GDAL 2.4.2, PROJ 5.2.0
library(maptools)
## Loading required package: sp
## Checking rgeos availability: TRUE
library(mapview)

base <-
  ATTdata$Station.Information %>% 
  st_as_sf(coords=c("Station.Longitude", "Station.Latitude"), crs = 4326) %>% 
  mapview(alpha.regions = 0, color = "grey", cex = 3, legend = F, homebutton = F, layer.name = "Receiver Stations")


## Now lets use the COA() function to calculate hourly centers of activity
COAdata <- 
    COA(ATTdata, 
        timestep = 60) ## timestep used to estimate centers of activity (in minutes)




## Now lets add our centers of activity for each individual on the map
COAmap <- 
  base + 
  COAdata %>% 
  st_as_sf(coords = c("Longitude.coa", "Latitude.coa"), crs = 4326) %>% 
  mapview(zcol = "Tag.ID", burst = T, alpha = 0, cex = 3, legend = F)

COAmap
proj<-CRS("+init=epsg:3174")
mcp_est <- HRSummary(COAdata,
                     projCRS = proj,
                     type = "MCP",
                     cont = 100,
                     sub = "%Y-%W")
## 
  |                                                                       
  |                                                                 |   0%
  |                                                                       
  |======================                                           |  33%
  |                                                                       
  |                                                                 |   0%
  |                                                                       
  |===========================================                      |  67%
  |                                                                       
  |                                                                 |   0%
  |                                                                       
  |=================================================================| 100%
head(mcp_est)
## $Overall
## # A tibble: 3 x 10
##   Tag.ID Sci.Name Common.Name Tag.Project Release.Date Tag.Life Sex   Bio  
##   <chr>  <fct>    <fct>       <fct>       <date>          <int> <fct> <fct>
## 1 22     Sander … walleye     HECWL       2012-03-27         NA FEMA… <NA> 
## 2 23     Sander … walleye     HECWL       2012-03-27         NA FEMA… <NA> 
## 3 153    Sander … walleye     HECWL       2012-03-20         NA FEMA… <NA> 
## # … with 2 more variables: Number.of.Detections <int>, MCP.100 <dbl>
## 
## $Subsetted
## # A tibble: 35 x 11
##    Tag.ID subset Sci.Name Common.Name Tag.Project Release.Date Tag.Life
##    <chr>  <chr>  <fct>    <fct>       <fct>       <date>          <int>
##  1 22     2012-… Sander … walleye     HECWL       2012-03-27         NA
##  2 22     2012-… Sander … walleye     HECWL       2012-03-27         NA
##  3 22     2012-… Sander … walleye     HECWL       2012-03-27         NA
##  4 22     2013-… Sander … walleye     HECWL       2012-03-27         NA
##  5 22     2013-… Sander … walleye     HECWL       2012-03-27         NA
##  6 22     2013-… Sander … walleye     HECWL       2012-03-27         NA
##  7 22     2013-… Sander … walleye     HECWL       2012-03-27         NA
##  8 22     2013-… Sander … walleye     HECWL       2012-03-27         NA
##  9 23     2012-… Sander … walleye     HECWL       2012-03-27         NA
## 10 23     2012-… Sander … walleye     HECWL       2012-03-27         NA
## # … with 25 more rows, and 4 more variables: Sex <fct>, Bio <fct>,
## #   Number.of.Detections <int>, MCP.100 <dbl>
## Lets have a look at how 100% MCP area changes over the tagging period
mcp_est$Subset %>% 
  mutate(Date = as.Date(paste(subset, "01", sep = "-"), "%Y-%W-%w")) %>% 
  ggplot(aes(x = Date, y = MCP.100/1e6, colour = factor(Tag.ID))) +
  geom_line() + geom_point() +
  labs(x = "Date", y = "Area of MCP (km2)", colour = "Tag ID") +
  theme_linedraw()

BBkud_est <- HRSummary(COAdata,
                     projCRS = proj,
                     h = 200,
                     type = "BBKUD",
                     cont = c(20,50,95),
                     storepoly = TRUE,
                     sub = "%Y-%W")
## 
  |                                                                       
  |                                                                 |   0%
  |                                                                       
  |======================                                           |  33%
  |                                                                       
  |                                                                 |   0%
  |                                                                       
  |===========================================                      |  67%
  |                                                                       
  |                                                                 |   0%
  |                                                                       
  |=================================================================| 100%
summary(BBkud_est)
##                 Length Class  Mode
## Overall         12     tbl_df list
## Subsetted       13     tbl_df list
## Spatial.Objects  3     -none- list
# Overall Brownian bridge UD estimates for full tagging period
head(BBkud_est$Overall)
## # A tibble: 3 x 12
##   Tag.ID Sci.Name Common.Name Tag.Project Release.Date Tag.Life Sex   Bio  
##   <chr>  <fct>    <fct>       <fct>       <date>          <int> <fct> <fct>
## 1 22     Sander … walleye     HECWL       2012-03-27         NA FEMA… <NA> 
## 2 23     Sander … walleye     HECWL       2012-03-27         NA FEMA… <NA> 
## 3 153    Sander … walleye     HECWL       2012-03-20         NA FEMA… <NA> 
## # … with 4 more variables: Number.of.Detections <int>, BBKUD.20 <dbl>,
## #   BBKUD.50 <dbl>, BBKUD.95 <dbl>
# Weekly Brownian bridge UD estimates for the three tagged walleye
head(BBkud_est$Subsetted)
## # A tibble: 6 x 13
##   Tag.ID subset Sci.Name Common.Name Tag.Project Release.Date Tag.Life
##   <chr>  <chr>  <fct>    <fct>       <fct>       <date>          <int>
## 1 22     2012-… Sander … walleye     HECWL       2012-03-27         NA
## 2 22     2012-… Sander … walleye     HECWL       2012-03-27         NA
## 3 22     2012-… Sander … walleye     HECWL       2012-03-27         NA
## 4 22     2013-… Sander … walleye     HECWL       2012-03-27         NA
## 5 22     2013-… Sander … walleye     HECWL       2012-03-27         NA
## 6 22     2013-… Sander … walleye     HECWL       2012-03-27         NA
## # … with 6 more variables: Sex <fct>, Bio <fct>,
## #   Number.of.Detections <int>, BBKUD.20 <dbl>, BBKUD.50 <dbl>,
## #   BBKUD.95 <dbl>
fullstack <-
  unlist(BBkud_est$Spatial.Objects)[grep("*_full", names(unlist(BBkud_est$Spatial.Objects)))]

names(fullstack) <-
  unlist(lapply(strsplit(names(fullstack), "[.]"), `[[`, 1))
library(raster, verbose = FALSE, warn.conflicts = FALSE)

## Lets plot the overall BBKUD for Tag.ID `153`
fulltag <- fullstack$`153`
values(fulltag)[values(fulltag) > 95] <- NA

## Raster plot
plot(fulltag, zlim = c(0, 100))

## Mapview plot
base + 
  mapview(fulltag)
library(leaflet, verbose = FALSE, warn.conflicts = FALSE)

## Full KUD for all tagged animals
fullmap <- COAmap@map

for (i in 1:length(fullstack)) {
  tempras<-disaggregate(fullstack[[i]], fact=3, method='bilinear')
  values(tempras)[values(tempras) >95] <-NA
  fullmap <- 
    fullmap %>% 
    addRasterImage(tempras, opacity = 0.7, group = names(fullstack)[i])
}

fullmap <- 
  fullmap %>%
  addLayersControl(
    baseGroups = names(fullstack),
    overlayGroups = "Receiver Stations",
    options = layersControlOptions(collapsed = FALSE)
  )